Avalanche dynamics, surface roughening and self-organized criticality - experiments 

on a 3 dimensional pile of rice. 
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We present a two-dimensional system which exhibits features of self-organized criticality. The 
avalanches which occur on the surface of a pile of rice are found to exhibit finite size scaling in their 
probability distribution. The critical exponents are r = 1.21(2) for the avalanche size distribution 
and D — 1.99(2) for the cut-off size. Furthermore the geometry of the avalanches is studied leading 
to a fractal dimension of the active sites of ds = 1.58(2). Using a set of scaling relations, we can 
calculate the roughness exponent a — D — ds = 0.41(3) and the dynamic exponent z — D(2 — t) = 
1.56(8). This result is compared with that obtained from a power spectrum analysis of the surface 
roughness, which yields a = 0.42(3) and z = 1.5(1) in excellent agreement with those obtained from 
the scaling relations. 
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I. INTRODUCTION 

The concept of self-organized criticality (SOC) 1 
presents a simple way of modelling slowly-driven out of 
equilibrium systems. The interesting natural systems 
thought to exhibit SOC, such as rain-fall 2 , earthquakes 3 , 
economic markets 4 , biological evolution 5 , or the brain 6 , 
are difficult to study in a controlled experiment, and more 
simple systems have to be found with which the predic- 
tions of SOC can be tested quantitatively. Among the 
first toy systems studied, were sand-piles, however it has 
been shown that due to the appearance of e.g. kinetic 
effects, real sand does not behave as an SOC system 7 . 
Another system, that has shown indications of SOC in 
a controlled environment is the dynamics of vortices in 
type-II superconductors. Here, power-law scaling has 
been observed in the past 8 , however not all authors have 
found avalanche scaling 9 . This is probably due to the fact 
that internal avalanches are not measured in the usual 
setup 10 . Recently Behnia et al. have measured inter- 
nal avalanches using arrays of Hall-probes, where they do 
indeed find power-law behavior of the avalanches. The 
ultimate hall-mark of SOC, however, is the observation 
of finite size-scaling, indicating a true critical dynamics. 
So far, there are only two experimental systems, a one 
dimensional (Id) pile of rice 12 or a Id pile of steel balls 
with a random arrangement at its bottom layer 13 , which 
have been shown to exhibit finite size scaling in their 
avalanche distributions. 

The critical exponents obtained from a finite size scal- 
ing analysis, would yield more information on a SOC 
system, as was discussed by Paczuski et al. 1A for models 
in a stationary critical state. A system evolving through 
avalanches, which are distributed according to a power- 
law, will also show roughening dynamics of its surface. 
This can be seen as an indication to the origin of the 
abundance of fractal or self-affine structures in nature. 
Roughening dynamics has been studied extensively in the 
past, including the study of interfaces in porous media 15 , 
the growth of bacterial colonies 16 , the slow combustion 17 



or the rupture 18 of paper . These systems are well charac- 
terized experimentally and in many cases, the interface 
dynamics can be modelled analytically by the Kardar- 
Parisi-Zhang (KPZ) equation 19 . In one spatial dimen- 
sion, the KPZ equation is solved exactly for the case 
of white noise, such that good comparison between ex- 
periment and theory is possible. However, e.g. in the 
burning of paper, there are avalanches observed in the 
propagation of the front, which are reminiscent of SOC 
dynamics 20 . This has recently been addressed by Alava 
et al. 21 , who proposed a mapping between SOC models 
and the KPZ equation similar to that of Paczuski and 
Boettcher 22 . With this mapping, they have been able 
to obtain SOC from the KPZ equation, with the proper 
kind of noise 21 . 

Moreover, roughening processes have also been found 
in the experimental SOC systems discussed above. The 
profile of a Id rice pile is a rough interface 23 , as is the 
front of magnetic flux penetrating a thin film supercon- 
ducting sample 24 . However, the universal scaling rela- 
tions derived by Paczuski et al. 14 have not yet been tested 
experimentally. But the fact that both roughening phe- 
nomena and avalanche dynamics can be observed in the 
same experimental system, such as a 3d rice-pile, means 
that the universal connections between SOC and surface 
roughening can be tested experimentally. 

Here we present a two-dimensional (2d) system, the 
surface of a 3d pile of rice, showing both roughening be- 
havior and avalanche dynamics. The avalanche dynam- 
ics is studied in terms of the avalanche size distribution 
for different sizes of the field of view, L. This allows 
the observation of finite size scaling and the determina- 
tion of the critical exponents D, describing the depen- 
dence of the cutoff scale on L, and r, the exponent of 
the avalanche size distribution. The roughening behav- 
ior is studied via the power spectrum of the surface 25 , 
in both space and time, resulting in a determination of 
the roughness and growth exponents respectively 26 . The 
connection between these two phenomena is shown from 
the derivation of scaling relations between the different 
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exponents, where we obtain excellent agreement using 
the experimentally determined values. 

The experimental setup for the growth of the rice-pile, 
as well as the reconstruction method used to determine 
the structure of the pile and the size of the avalanches is 
discussed in section II. The avalanche size distributions, 
and their finite size scaling is presented in section IIIA, 
together with the determination of the critical exponents. 
In section IIIB, the surface roughness is analyzed and the 
necessary techniques are briefly introduced. Still in sec- 
tion IIIB, the scaling relations between the roughness ex- 
ponents and the critical exponents are introduced. These 
results are also put into a wider context and compared 
with results from KPZ roughening systems 27 . 



II. EXPERIMENTAL DETAILS 

The experiments were carried out on long grained rice 
with dimensions of typically ^2x2x7 mm 3 , similar to rice 
A of ref. 12 . The pile was grown from a uniform line source 
1 m wide. This uniform distribution was achieved via 
a custom built mechanical distributor based on a nail- 
board producing a binary distribution 28 . The actual 
setup consists of a board with an arrangement of tri- 
angles as shown in Fig. 1. This means that the possible 
pathways of the rice are continuously split at each level, 
such that at the bottom we end up with a row of 64 uni- 
formly distributed compartments. The uniformity of the 
distribution was measured to be of the order of 5%. 




FIG. 1. A schematic image of the setup. The distribution 
board can be seen on top, where rice is dropped from a sin- 
gle point and subsequently divided into even compartments. 
Within the wooden box bounding the rice pile, a reconstruc- 
tion of its surface is shown, as it is used in the further analysis. 

At the bottom of the distributor, the grains are slowed 
down by a sheet of plastic before they hit the top of the 
pile. The rice is fed to the top of the distribution board, 
from a point source, which drops rice at an average rate 
of ^5g/s. This corresponds to 1500 grains per image, 



which means that over the length of the line there are 1 
to 2 grains dropped per place. This rate is uniform over 
the time scale between two images. 

Once a rice pile is grown, we measure the further evo- 
lution of the surface coordinates using a specially devel- 
oped real-time technique based on the projection of a set 
of lines. In order to increase the spatial resolution in 
the direction of growth, we use a set of lines in the base 
colors (red, green and blue), which can be easily sepa- 
rated digitally for the analysis. The lines are projected 
at a right angle to the average surface of the pile and 
observed at 45° to both the direction of projection as 
well as the surface. This is illustrated in Fig. 2, where 
a raw image is shown together with the lines extracted 
from it. From the distortions of the lines observed in this 
way, the coordinates in 3d space of the pile surface can 
be calculated by simple geometry 28 . With this recon- 
struction technique, the 3d coordinates of the whole field 
of view can be determined with an accuracy of 2-3 mm, 
as we have checked on a number of test surfaces. Using 
both structured and smooth surfaces, the resolution and 
accuracy were determined independently. Since the reso- 
lution is roughly matched with the size of the grains, the 
method is well suited for the present purpose. 




FIG. 2. A raw image of the rice pile with the identified 
lines. From the distortions of the lines, the surface geometry 
is obtained. The dark spot in the image indicates a reference 
point. 
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In a single experimental run, the growth is studied for 
a period of ~4 hours, with a picture taken every 30 s. 
Thus an experiment consists of ~ 480 images. The pic- 
tures are taken with a high resolution digital camera, 
having a resolution of 2048x1536 pixels. For every time 
step, the surface structure is reconstructed, which gives 
information about the roughening process, as well as the 
avalanche dynamics. 
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FIG. 3. A snapshot of an avalanche, given by the height dif- 
ference of two subsequent images. In the picture, the height 
is given by the grey-scale with white corresponding to A/i = 
and black to A/i = 15 mm. Integrating the height difference 
Ah(x, y) over the area gives the size of the avalanche in terms 
of a volume, or in terms of the number of moved grains from 
the volume occupied by a grain V gr ain = 35 mm 3 . 

The size and shape of the avalanches can be deter- 
mined from the height difference of the surface between 
two consecutive images. The overall growth of the pile 
is subtracted from this difference, however this correc- 
tion is negligible and does not influence the results. This 
is shown in Fig. 3, where the height difference A/i(cc, y) 
is shown for a medium size avalanche. This allows the 
study of internal avalanches instead of just the off-edge 
ones, which has proven to be important in previous stud- 
ies (for instance on superconducting vortex avalanches 11 
and Id piles of rice 12 ). In order to obtain the size of 
such an avalanche, \ Ah\ is integrated over the area, which 
yields the displaced volume, AV, corresponding to the 
size of the avalanche in mm 3 . In order to use natural 
units, we will in the following measure the size of the 
avalanche in terms of the number of toppled grains, which 
is s = AV /Vg ra i n , where Vgrain — 35 mm 3 is the volume 
of a single grain of rice. 

The data discussed in this article was obtained in three 
separate experiments, with a total of 1330 images. 



III. RESULTS 

In this section the experimental results are presented. 
We first show that there is finite size scaling in the 



avalanche statistics, which indicates the appearance of 
a critical state in the system. Second, we characterize 
the surface roughness of the rice pile in 2d, both in space 
and time. The two characteristic exponents, a for the 
spatial behavior, and z for the dynamics 26 can also be 
obtained from scaling relations of the critical exponents 
of the avalanche statistics and are in excellent agreement 
with those determined from a standard power-spectrum 
analysis. 



A. Avalanche Dynamics 

The time dependence of the displaced volume, AV, in 
each time step is shown in Fig. 4, for all three experi- 
ments discussed here. In the figure, the data are shown 
for avalanche sizes integrated over the total field of view 
of 600x600 mm 2 . 
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FIG. 4. The time evolution of the integrated height dif- 
ference for the three experimental runs studied. As can be 
seen, times of relative rest are punctuated by large avalanches, 
which appear intermittently. When studying the distribution 
of avalanche sizes (see Fig. 5), it is found that there is no 
intrinsic size to the avalanches, but that they are distributed 
according to a power-law. The different experimental runs 
are indicated by arrows. 

A histogram of these data, giving the avalanche size 
distribution is shown by the diamonds in Fig. 5a). In 
that figure, size distributions of subsets of the data cor- 
responding to smaller fields of view (L = 50, 150 and 300 
mm) are also shown. Due to experimental resolution, the 
smallest avalanches measured for each subset depend on 
the size of that subset. The range of sizes from 50x50 
mm 2 to 600x600 mm 2 spans more than a decade and all 
of the data taken together show a power-law scaling of 
the size distribution over three decades, with an exponent 
of T = 1.20(5). 

These data can however be scaled to fall onto a single 
curve, as shown in Fig 5b). Here the avalanche sizes are 
scaled with L~ D and the probabilities are scaled with s T . 
The good curve collapse visible in Fig 5b) indicates the 
presence of finite size scaling in the data, one of the hall- 
marks of critical behavior in a system and thus evidence 



of the appearance of SOC in the 3d rice pile. 




FIG. 5. a) The unsealed size distribution functions for the 
avalanches from different subsets of the experiments corre- 
sponding to sizes of L = 50, 150, 300 and 600 mm. Given 
the asymptotic dependence of all them together, there is 
power-law behavior over at least three decades, with an ex- 
ponent of r = 1.20(5). b) The same data, scaled to produce 
a curve collapse. The sizes of the avalanches are scaled with 
L~ D and the probabilities are scaled with s T . The values 
used to obtain the best curve collapse are r = 1.21(2) and D 
= 1.99(2). 

Due to the finite size scaling found from the curve col- 
lapse, the avalanche size distribution as a function of sys- 
tem size can be written as a function of one parameter 
only: 

P( s ,L)=s-Tf(JL 5 ), (1) 

where f(x) is constant up to a cutoff scale 

s co ^L D . (2) 

The exponents r and D used to obtain the curve collapse 
in Fig. 5b) were r = 1.21(2) and D = 1.99(2). We note 
here that in usual finite size scaling, separate experiments 
are used instead of subsets. However, as we have tested 
on simulations of a 2d version of the Oslo- model 29 , finite 
size scaling using subsets yields the same exponents, but 



is much more easily implemented experimentally. As an 
additional test, the avalanche dimension can be deter- 
mined directly using a box-counting method in 3d. This 
yields D — 2.05(10) consistent with the result from finite 
size scaling. 

B. Surface Roughening 

The roughness of a surface can be characterized in dif- 
ferent ways. Most commonly, this is done via the width 
of the interface as given by its standard deviation 26 . For 
a 2d surface, this is given by 

™ 2 (L) = ±jr(h(i,j)-{h)f, (3) 

where (h) is the average profile of the surface height. 
When the roughness of the surface has reached satura- 
tion, its value w sa t, will scale with the system size as a 
power-law, w sa t oc L a , where a is the roughness expo- 
nent. 
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FIG. 6. Determination of the roughness exponent from a 
power-spectrum analysis as well as from the correlation func- 
tion C(x) (squares). After subtracting the average slope of 
the pile, its 2d power spectrum is calculated. After taking a 
radial average, it is integrated over k-space in order to give 
the distribution function a, which follows a power-law de- 
pendence a oc k~ a , with a roughness exponent a = 0.42(3). 
Similarly, the correlation function C(x) oc x a (squares) is cal- 
culated from the pile surface, resulting in the same exponent 
as that determined from the distribution function. 

A similar type of analysis can be achieved via the power 
spectrum of the surface. Here, one has the advantage 
that better statistics can be obtained, since the whole 
surface can be used at once 25 . We determined the power 
spectrum from a radial average of the 2d Fourier trans- 
form of the surface 

S(k) = \h(k X ,ky)\ 2 , (4) 

where k = (k x + ky) 1 / 2 . The Fourier transforms where 
performed via an FFT algorithm, which is why a data 



subset corresponding to a power of two (512x512) was 
used. From the power spectrum we then determine the 
distribution function a(k) 



D-dt 



(7) 



S{K)Kdn, 



(5) 



from which it can be shown 25 that cr(27r/fc) = w(£), where 
I = 2ir/k is the length scale over which the width is 
calculated. The roughness exponent a can therefore be 
reliably obtained from a determination of a(k). An al- 
ternative way of obtaining the roughness exponent is via 
the two point correlation function 26 . 



C(x, t) = ((h(£, t) -h(x + £,t + t)Y) 1t 



1/2 



(6) 



the radial average of which obeys the same scaling be- 
havior as the distribution function 26 . 

In Fig. 6, the correlation and distribution functions 
for the rice-pile surface is shown, averaged over all time 
steps of all experiments (1330 images). The part of the 
surface studied consists of 512x512 mm in the center of 
the surface laterally (along the horizontal direction) and 
starting from a height of 200 mm. As can be seen from 
the double logarithmic plot, there is power-law scaling 
with an exponent, a = 0.42(3). 
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FIG. 7. Determination of the fractal dimension of the 
avalanches. Over all experiments, avalanches with a size two 
standard deviations bigger than the average were studied. 
The contours are derived from a thresholding at the value of 
the mean height difference for each image. The fractal dimen- 
sion is subsequently determined from a simple box counting 
method. 

The scaling behavior of the interface width can, how- 
ever, also be determined from the cut-off size of the 
avalanches in a system with size L. The size of such 
an avalanche will be proportional to the maximal height 
difference, given by the saturation width w sa t, times the 
maximal area of the avalanche. As can be seen from 
Fig. 3, the avalanche-shape is a fractal, which means that 
s co oc w sa tL dn , where ds is the fractal dimension of the 
cluster of displaced grains. Equating this expression with 
Eq. 2, one obtains w sa t oc L D ~ ds , which implies the scal- 
ing relation 



Such a scaling relation between the roughness exponent 
and the exponents characterizing the avalanches was also 
derived by Paczuski et al. 14 from different arguments. 

In order to determine ds in an accurate manner, we 
studied avalanches more than two standard deviations 
bigger than the average size. This corresponds to a sub- 
set of ~ 100 images from all experiments. After thresh- 
olding the height difference fields Ah(x,y) at their mean 
value, we applied a simple box counting method 30 to the 
resulting clusters. This is shown in Fig. 7, where the 
number active pixels in a box of given size is shown as a 
function of the box size. The result is a power law with 
an exponent of -1.58(2), indicating a fractal dimension of 
c?b = 1.58(2). Inserting this value into scaling relation 
(7), the roughness exponent is found to be a = 0.41(3), 
in excellent agreement with that determined from the 
roughness analysis. 
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FIG. 8. Determination of the growth exponent from a 
power-spectrum analysis. Here, the distribution function 
a has been determined from the time dependence of the 
height of each pixel. The data show a power-law dependence 
a oc uj' 13 , with a growth exponent (3 = 0.28(3). The same re- 
sult is obtained via the correlation function, indicated by the 
open symbols. From the growth and roughness exponents, 
the dynamic exponent z = a//3 = 1.5(1) can be determined. 

The dynamics of the roughening process is likewise an- 
alyzed via the distribution function. In contrast to the 
roughness analysis above, we now determine the distribu- 
tion function for the time dependence hij(t) — (hij)t for 
each pixel, where (-) t denotes the average over the dura- 
tion of the experiment. Again, the distribution function, 
given by 



a 2 (uj) = / \h{uj)\ 2 dw, 



(8) 



is equal to the momentary width, <T(27r/w) = w(t). Thus 
the growth exponent (3, describing the scaling of the 
width with time, can be determined from the distribu- 
tion function. Again, the correlation function shows the 



same scaling behavior in time, allowing a separate deter- 
mination of (3. Both results are shown in Fig. 8, where it 
can be seen that there is good power-law scaling of the 
distribution function, as well as the correlation function, 
over two decades with an exponent of /3 = 0.28(3). An- 
other way of describing the dynamics of an interface is 
via the dynamic exponent z, which describes the scaling 
of the saturation time with the system size 



t x oc L z 



(9) 



It can be easily obtained 26 from the roughness and 
growth exponents via 



(10) 



From the values determined above we thus obtain z = 
1.5(1). 

The dynamic exponent, z, can also be derived from 
a scaling relation using the critical exponents of the 
avalanche dynamics. The saturation time, t x , will be 
roughly given by the time it takes for an avalanche of 
the cut-off size, s co , to appear. Since material is added 
to the system at a constant rate, the number of grains 
added until a cut-off avalanche occurs is proportional to 
the cross-over time t x . On the other hand, the material 
added not only result in an avalanche of size s co , but will 
also be lost by smaller avalanches, such that the total 
material necessary to obtain an avalanche of size s co can 
be estimated from: 



sP(s)ds oc s 2 



(11) 



From Eq. 2, we obtain for the cross-over time t x oc 
£D(2-t) with Eq. 9, this immediately leads to the scal- 
ing relation 



D(2 



(12) 



which was also derived by Paczuski et al. from more 
general arguments. Using the values for D and r de- 
termined above from the curve collapse of the avalanche 
size distributions, we obtain z = 1.56(8), which is again 
in very good agreement with the result of the power spec- 
trum analysis, z — 1.5(1). Note also that the exponents 
a and z fulfill the KPZ scaling relation 19 a + z = 2, which 
is valid independent of the dimension of the system and 
depends only on the fact that the growth is driven by 
height gradients. Furthermore, the roughness exponent 
a = 0.42(3) is in good agreement with a numerical de- 
termination of the behavior of the 2d KPZ equation 27 . 
The connection between roughening and SOC will be dis- 
cussed in more detail below. 



IV. CONCLUSIONS 

In conclusion, we have presented results on both the 
avalanche and the roughening behavior of an experimen- 
tal SOC system. In addition, we presented simple argu- 
ments for universal scaling relations, derived by Paczuski 



et al. 1A on general grounds, connecting the avalanche and 
roughening behavior. We obtain quantitative agreement 
of the experimental exponents characterizing the rough- 
ness and those describing the avalanche statistics. This 
means that the roughening of the surface of the pile is 
governed by its underlying avalanche dynamics, which 
was already conjectured from models of interface depin- 
ning. Earlier experimental studies 12 ' 23 have considered 
the avalanche dynamics and the roughening behavior of 
the Id version of the present system, where finite size 
scaling was also found 12 . The above scaling relations, 
showing quantitatively the connection between the two 
phenomena were not previously tested. 

The fact that KPZ dynamics is observed in a SOC 
system, as indicated by the fact that the roughness and 
dynamic exponents obtained here fulfill the KPZ scaling 
relation, is surprising. This is mostly because in SOC 
systems, the critical state is self organized, whereas in 
KPZ systems, the roughening is put into the dynamics 
by obeying the proper symmetries. However, it has been 
argued by Alava et al. 21 , that there is a mapping of SOC 
models, which naturally lead to the necessary symme- 
tries. In this context, it should also be noted that exper- 
iments on KPZ systems such as burning paper have also 
observed 'avalanches' in the advance of fronts 20 . 
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